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Abstract 

We consider stagnation point flow away from a wall for creeping flow of dilute poly- 
mer solutions. For a simplified flow geometry, we explicitly show that a narrow 
region of strong polymer extension (a birefringent strand) forms downstream of the 
stagnation point in the UCM model and extensions, like the FENE-P model. These 
strands are associated with the existence of an essential singularity in the stresses, 
which is induced by the fact that the stagnation point makes the convective term 
in the constitutive equation into a singular point. We argue that the mechanism 
is quite general, so that all flows that have a separatrix going away from the stag- 
nation point exhibit some singular behaviour. These findings are the counterpart 
for wall stagnation points of the recently discovered singular behaviour in purely 
elongational flows: the underlying mechanism is the same while the different nature 
of the singular stress behaviour reflects the different form of the velocity expansion 
close to the stagnation point. 
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1 Introduction 



Extensional flows of polymer solutions and melts occur in many industrial 
polymer processing operations, and hence such flows have been studied for 
decades [TH2] . Recently however, interest in extensional flows was renewed 
by observations of steady and unstable continuous flow in microfluidic de- 
vices [3]Hf5] , and it was realized only recently that extensional flows are prone 
to the formation of singularities and non-analytic structures in the stress fields 
[oTHIg] . Depending on the Deborah number and the model used, these stress 
singularities may take various forms. For purely extensional flow in contin- 
uum models that describe infinitely extensible polymer chains (such as the 
upper convected Maxwell model (UCM) and the Oldroyd-B model [Tf2"1l9] ) 
the stresses can have power law spatial behaviour with a finite limit at the 
centre line, or they even have power law divergencies. For models that are 
based on finitely extensible chains, divergencies are cut off at some scale, but 
singular behaviour of the stress gradients may persist [6.7.8J. Such singular 
behaviour may have important implications for numerical simulations of ex- 
tensional flows, since it leads to structures with a very small length scale. 
Indeed it is known that for many such flows, numerical schemes break down 
at only moderate flow rates (Deborah numbers of order unity). 

The question quite naturally comes up whether singular behaviour near special 
points is the rule rather than the exception. We argue in this Communication 
that the latter is the case and demonstrate this for a simplified case where all 
calculations can be done analytically, so that the emergence of the singular 
behaviour can be followed explicitly. 

The reason to expect singular behaviour near special points where the velocity 
vanishes — even though the geometry is not sing - is actually very 

simple. For steady flow, the only derivative terms of the stress T in UCM- 
type constitutive equations come from the convective term (v • V)T. The 
points where v vanishes — the stagnation point in elongational flow or in the 
wall stagnation point flow considered here — thus translate into a singular 
point [10] of the partial differential equation obtained from the constitutive 
equation for the stress. Close to the singular point, the lowest order terms 
in the expansion of v are often fixed by simple symmetry considerations and 
boundary conditions, if applicable. So the nature of the dominant singularity 
at the singular point is generally fixed independent of the precise details of the 
model. Further away from the singularity, the behaviour will typically depend 
on the details of the flow profile. All these features are well illustrated by the 
analysis below. As stated, we focus on a simple case where the calculations 



Of course, at sharp corners where the flow field itself is singular, this singular 
behaviour carries over to the stresses. 
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Fig. 1. Stagnation flow (a) in a wake, (b) approximated by a flow near a flat wall. In 
(b) the formation of a birefringent strand is qualitatively indicated by the shaded 
area. 

can all be done analytically, but the scenario holds generally for complex more 
realistic flows and we suspect this mechanism of advection to be at the origin 
of the formation of birefringent strands. 

We focus on wall stagnation point flows where the flow is away from the 
wall. Examples of this are flows in the wake of a falling sphere or of a fixed 
cylinder, as shown schematically in Fig. [1] (a). In particular, the flow past 
a fixed cylinder or sphere in a channel has become a benchmark problem 
for numerical modelling of viscoelastic constitutive equations [TTf 121113] . It is 
known that in such flows a narrow region of high polymer extension may form, 
a so-called birefringent strand [T^IT5] . This region starts at a small but finite 
distance downstream from the stagnation point, as indicated schematically in 
Fig. [1] (b), where a flow near a flat wall is depicted. 

In this work, we consider a strictly two-dimensional version of this flow, with a 
simplified, fixed velocity field obeying the basic symmetry of a stagnation point 
at a wall (cf. [H]). We analyse this case in detail for the UCM model [H21I9] 
but also discuss in the end the qualitative changes that occur for a FENE-P 
model. 

Unlike the case of steady purely extensional flow, which was analysed previ- 
ously [6f7f8] , the extension of the polymers does not diverge for any extension 
rate. We find that a thin birefringent strand forms, with a singularity at its 
centre. As argued above, notwithstanding the simplifications we make in ob- 
taining this result, we believe that the analysis makes it clear how singular 
behaviour emerges in general. 

Unfortunately, our results cannot immediately be compared quantitatively 
with experiments or numerical computations on realistic cases like flow past 
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a cylinder [T2lfT3] . First, one should keep in mind that in such situations, 
there may be two sources of (near) singular behaviour: besides the one we 
analysed here, dominated by the symmetry and boundary conditions of the 
velocity field near the wall stagnation point, in viscoelastic flow past a cylinder 
large stress fields are already built up at the sides of the cylinder, where the 
flow is mostly along the cylinder. These shear stresses are advected toward 
the rear stagnation point. This effect is clearly not present in the simplified 
geometry that we consider. Second, our analysis is based on taking a fixed 
velocity field obeying the basic symmetry, and we show how this leads to an 
essential singularity in the stress. In reality, of course, the velocity and stresses 
are coupled indirectly through the momentum balance equation. Through this 
coupling, the velocity field will also be affected in the region of large stress 
gradients near the stress singularity. Since the symmetry and expansion of the 
velocity near the stagnation point cannot change in lowest order, we expect 
that there is an intermediate flow regime where the basic structure of the 
singularity is not changed dramatically. This assumption is further supported 
by recent simulations of a two-dimensional cross-slot flow by Poole et al. [5]. 
There, the velocity profiles remained smooth even when the flow changed its 
symmetry (the new type of purely elastic instability discovered by Arratia 
et al. [1]), while stresses exhibit the typical singular structure similar to the 
one discussed in [6f8] . At the same time, numerical studies suggest that at 
sufficiently large flow rates, this nonlinear coupling can become so strong that 
there may be no steady state flow solution past a cylinder for Deborah numbers 
of order unity [T^fTB] . The coupling and this effect are, unfortunately, beyond 
the present approximation. 



The layout of this paper is as follows. In section 2 we introduce the flow 
geometry and the models, and we briefly recapture similarity solutions for 
UCM found by other authors [TTyi~8] . We calculate analogous solutions for a 
simplified version of this flow, where we fix the velocity field, for UCM and 
FENE-P. In section 3 we consider more realistic boundary conditions, and 
we solve the constitutive equations analytically. In section 4 we consider the 
resulting stress field (extension field) in more detail, showing that we find a 
narrow region of high polymer extension, with a non-analytic stress profile at 
the centre of the strand. We then discuss these results in the light of more 
realistic flow profiles, and we conclude by discussing the relevance of these 
results for computational and experimental work. 
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2 Simplified stagnation flow of a UCM fluid 



We consider incompressible planar stagnation flow of a UCM fluid without 
inertia (creeping flow). The UCM constitutive equation for steady flow is [T] 



T + A (v • V)T - (Vv) J T - T (Vv) = 7? Vv + (Vv 



where A is the relaxation time of the fluid and 77 is the Newtonian viscosity. 
The momentum balance for creeping flow is 



V • T — Vp = 0, 

where p is the pressure. Incompressibility is given by 

V • v = 0. 



(2) 



(3) 



The planar stagnation flow geometry we consider is depicted in Fig. [2J We take 
the vertical direction as bounded, with length £. Because of the solid wall, the 
boundary condition at the wall (y = 0) is v = 0. At y — £, we impose v y = V, 
with V > 0. For the velocity field, a similarity solution then exists, which is 
of the form [17118] 

v x = -xip'(y) and v y = ij)(y), (4) 

where the boundary conditions imply i/)(Q) = and ip'(0) = 0. Inspired by this 
solution, we take a fixed velocity field that satisfies the boundary conditions 
and that would correspond to the lowest-order approximation near the wall: 

v x = -2(V/l 2 )xy and v y = (V/£ 2 )y 2 . (5) 

Note that these terms are the lowest order analytic terms in an expansion in 
x and y away from a symmetric stagnation point at the wall. This is why our 
analysis illustrates the emergence of singular behaviour in more general cases 
as well. 

Let us now insert this velocity field into the constitutive equation. The result- 
ing stress field will not in general satisfy momentum balance, but it does yield 
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Fig. 2. Planar stagnation flow. 
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a valid solution for Newtonian creeping flow. This is sometimes referred to as 
"Newtonian kinematics" ; in the Oldroyd-B extension of the UCM model, this 
would be a reasonable approximation in the dilute limit, (5 -C 1, in which the 
polymer stresses do not influence the flow [rSlfTB] . 

We can rescale the quantities appearing in the equations: length is scaled with 
£, velocity with V, and stress with rjV/£. The velocity field becomes 



-2xy and v y 



y 



(6) 



with < y < 1, and for the constitutive equation we obtain a dimensionless 
form 



T + De [(v • V)T - (Vv) T T - T (Vv)] = Vv + (Vv) T (7) 

Here we introduced the Deborah numbeJ" 2 ! 

De = XV/l (8) 

We can then insert the rescaled velocity into the constitutive equation and 
solve for the stresses. We obtain equations for the components of the stress 
tensor, and we observe that the equation for the yy component decouples from 
the other equations, because dv y /dx is identically zero for this velocity field. 
The equation for the yy component becomes 



T 



illl 



De 



8T 
-2xy a -j^ 

ox 



+ y 



2 u± yy 
dy 



4yT, 



yy 



Ay. 



(9) 



Let us first analyse the solution of this equation under the often-made assump- 
tion that T yy is constant in x [17J and then analyse why this solution misses 
an important part of the physics. One might naively think that this solution 
can be seen as the first term in a series expansion in powers of x [18J, but 
as we shall see this assumption itself is incorrect: singular terms are typically 
generated by the stagnation point flow. Defining 

T yy (x,y) = Y(y), (10) 

we solve 

Y(y) + De [y 2 Y' (y) - AyY(y)] = Ay. (11) 
This equation allows an exact solution, and we find 

Y gen (y) =Ay + 12Dey 2 + 24DeV + 24DeV + Cy 4 exp[l/(Dey)], (12) 



2 One should keep in mind that our Deborah number cannot directly be compared 
with the one used in studies of flow past a cylinder. 
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Fig. 3. Planar stagnation flow on the finite rectangular domain analysed in this 
paper. 



where the last term corresponds to the homogeneous solution of the differen- 
tial equation, and C is a constant of integration. This term clearly leads to 
unphysical results, as it implies — > oo as y — > + . Hence, if one thinks 
of Y(y) as the first term in a regular series expansion in x, we should discard 
it and keep only the particular solution 

Y (y) = Ay + 12De?/ 2 + 24DeV + 24DeV. (13) 



A few remarks about this solution are in order: (i) It does not in general 
satisfy momentum balance, but it is qualitatively similar to the more accurate 
similarity solution of Oztekin et al. [IB]- (ii) By forcing the solution to be 
independent of x, we find a solution that can only be seen as an approximation 
that is valid in a small range around x — 0. (Hi) We assume the same functional 
form of the velocity field for all De. (iv) The solution does not diverge at finite 
y for any De. This is in contrast with purely extensional flow, where stresses 
diverge for De > 1/2. 

Let us now discuss the shortcoming of the above line of analysis. The solution 
is an admissible solution if we work in the infinite domain x — > ±oo. How- 
ever, if we work in a finite domain, with boundary conditions at some finite 
x = ±L, say, then the solution is only relevant if the stresses on this boundary 
are precisely consistent with fflBl . In general, the boundary stresses are incom- 
patible with this expression and the flow drawn in Fig. [2]advects all deviations 
from this expression towards the stagnation point. While in the above analysis 
it appears as if the analysis close to the stagnation point dictates what the 
proper boundary condition far to the left and right should be, in reality we 
have to analyse what happens with the stresses as they are advected from the 
left and right to the stagnation point. That dictates the behaviour there, not 
the other way around! The analysis below does show that singular behaviour 
is picked up through this advection from the boundaries, as one might expect. 
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3 Inflow boundary conditions and an explicit solution 



To back up the above observations, we now allow explicit x-dependence of 
the stress field without including higher orders in x in the velocity field. We 
can then match more realistic inflow boundary conditions by including ho- 
mogeneous parts of the solution, and we shall see that this causes qualitative 
changes in the stress field. 

To do so, we modify the geometry by explicitly taking a finite domain in the 
x direction, with — 1 < x < 1, see Fig. [3j We observe that the equation is 
purely advective. This implies that is it a first-order differential equation and 
that there is no "interaction" across streamlines. We can therefore split the 
domain in two separate parts, one for x < and one for x > 0, with the line 
x = acting as the separatrix between the two domains. As in the case of 
pure extensional flow, this allows for non-analyticity on the line x = 0. [6117118] . 
On each of the two subdomains, we can impose independent inflow boundary 
conditions. We will restrict ourselves to the domain < x < 1, and we will 
assume that the negative domain is the mirror image of this domain. On the 
lateral boundary x = +1 we impose the boundary condition T yy (x — l,y) — 
Y h {y) with Y b (0) = 0. 

We thus consider again the differential equation (j5J): 



Y(x,y) + De 



ox ay 



Ay. (14) 



In the language of local analysis of differential equations, this partial differen- 
tial equation has an irregular singular point at x = y = [TO] . 

It is easy to check that the partial differential terms in the equation together 
give identically zero if Y is a function of xy 2 only. In other words, the form 

Y(x,y) = f(xy 2 ), (15) 

is a zero mode of the differential operator since 

dfjxy 2 ) 2 df(xy 2 ) 

-2xy — + y — = (16) 

ox oy 

for any function /. 

Now, due to the linearity of the equation, the solution of the equation is 
the sum of the particular solution Y given in (1131) plus an arbitrary solution 
Y horn (x,y) of the homogeneous equation, that is 

Y(x,y) = Y (y)+Y ham (x,y). (17) 
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Moreover, if we assume Yhom to be of the form 

Y hom (x,y) = g(y)h(xy 2 ), (18) 

then inserting this into the homogeneous part of the equation (with right hand 
side equal to zero) effectively gives an equation for g(y) that is identical to 
the homogeneous part of Eq. ffTTj) . For g we thus recover the homogeneous 
solution of ffTTT) . that is, the last term of (fT2"l) . Thus, we have 

Y(x,y) = Y (y)+y 4 exp[l/(Dey)]h(xy 2 ) (19) 



We can now impose the boundary condition Y(x = l,y) = Y^(y) at x = 1: 
requiring that Y(x = l,y) given by f|T9|) equals Y^(y), immediately gives 

h(y 2 ) = Y *M- Y °M eX p[-l/(De 2/ )]. (20) 

Since h(xy 2 ) is a function of the product xy 2 only, we now know h for all x 
and y. This finally gives the general full solution 



F(x,y) = Y {y) + — * '— * '- exp 



1 ~ 1/yg 
Dey 



(21) 



Note that the structure is precisely as we envisioned: the deviations of the 
x-independent solution Yo(y) from the stress boundary condition Y^(y) are 
convected towards the stagnation point, and the stresses there are largely 
dominated by the singularity. 

For the flow that we consider, we use the procedure above to match the stress 
field to the inflow boundary condition. For simplicity we take this to be 

Yb(v) = 0, (22) 

but since the behaviour for small x is dominated by the singular term, other 
choices lead to similar conclusions. For De = 1 we then find the stress field in 
Fig. HI Since we consider the ultradilute limit, one may also think of this as 
the polymer extension in the yy direction (see below). 

We clearly see that a narrow region of high extensional stress (polymer exten- 
sion) forms in the centre of the domain. The solution we find does not diverge 
for any De, nor does any gradient diverge — in fact, the x-derivatives are zero 
to any order on the centre line, and we have an essential singularity there. 
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(b) 



(c) 



Fig. 4. Polymer extension in a simplified wall stagnation flow for ultradilute UCM 
fluid, for De = 1.0: (a) three-dimensional plot of the extension as a function of x and 
y, (b) contour plot of the extension, showing the shape of the birefringent strand, 
(c) cross-section of the extension profile for y = 0.25,0.5,0.75, 1.0 






Fig. 5. Polymer extension in a simplified wall stagnation flow for ultradilute FENE-P 
fluid, for De = 5.0, and L 2 = 10.0, 100.0, 1000.0. 



4 Distance from the wall 



The stress profile Y (y), Eq. (fT3|) . for a UCM fluid has no intrinsic length scale, 
as the UCM model has no characteristic stress or extension. Therefore, we can- 
not obtain a meaningful "distance from the wall" in that model. Models with 
finite extensibility, such as the well-known FENE-P model, are required to do 
this. In previous work we argued [8] that for extensional flow a fair approxima- 
tion to the FENE-P model can be obtained by restricting to the yy extensional 
component of the stress tensor, and making the simple approximation that the 
polymer stretching follows the UCM rheology up to a certain maximum ex- 
tensional stress after which the stress does not increase any further (this is 
sometimes referred to as a linear-locked approximation [19 20J). This maximal 
stress should correspond to the maximal extension L in the FENE-P model, 
for which the trace of the conformation tensor A is L 2 . The UCM stress tensor 
is related to the FENE-P conformation tensor as A = 1 + De T (this can be 
seen by taking L — > oo in the FENE-P constitutive equations, and comparing 
this to the UCM constitutive equations). Restricting to the yy component, 
we find 1 + DeT™ ax = L 2 . In Fig. Owe show stress profiles for the simplified 
velocity field, for the UCM model and a FENE-P-like model that is restricted 
to the yy component. 

In this approximation and for a fixed velocity field, we can thus estimate the 
distance of the strand from the wall, as well as the width of the strand, by 
considering the locus of points where the UCM extension 1 + De T yy reaches 
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L 2 . For our fixed velocity field, it is straightforward to find scaling relations 
for the distance and the width as functions of De and L 2 . From Eq. (|T3|) we 
see that the distance y^ from the wall is given in our approximation by 

1 + De y + 12 De 2 yl + 24 De 3 yl + 24 De 4 = L 2 . (23) 

This means that yo scales as 1/ De for fixed L 2 ^> 1 and all De. For fixed De 
we find that y should scale as (L 2 ) 1//4 = \f~L for L 2 — > oo. We conclude that 

yo = ^^rrJ(L 2 ) with f(L 2 ) < 1 and f(L 2 - oo) = 1. (24) 

This asymptotic result cannot be easily extended to other, more realistic ve- 
locity fields. Note that in the limit of large L 2 , the distance yo will become 
larger than 1, and it falls outside the domain that we assumed, < y < 1. The 
general conclusion that the distance increases with L 2 and decreases with De 
is of course physically reasonable: for larger L 2 at fixed De it takes a longer 
distance to fully stretch the polymers, while for larger De at given L 2 , the 
stretching occurs more rapidly and one may expect that the distance required 
for full stretching decreases, at least on the centre line of the flow. 



5 Width of the strand 

The width of the strand is a more subtle issue. Since the stress is bounded, 
and therefore also the UCM extension, we can define a meaningful length 
scale even in the UCM model, as can be seen from the left panel of Fig. For 
example, we can take the full width at half maximum (FWHM) of the peak. 
In that case, the width decreases as a function of De. We can also look at our 
FENE-P-like approximation and use the contour line where the profile reaches 
a given L 2 . A numerical analysis shows that the FWHM decreases as 1 / De , 
while the approximate FENE-P width increases monotonically and saturates, 
see the right panel of Fig. [6j 

These asymptotic results can be explained from the exact solution, Eq. fl2~T|) . 
The FWHM is mostly determined by the behaviour close to the centre of the 
strand, where the factor exp[— l/(De \/x)} dominates: as x — > it vanishes 
more rapidly than the factors in front of it diverge. The typical length scale 
x for this exponential factor is given by -^/xoDe ~ 1, or x oc 1/De 2 , and 
thus, the FWHM also scales as 1/De 2 . For the FENE-P-like width we make 
the observation that the UCM extension as a function of x actually converges 
pointwise as De — > oo: 

lim (1 + DeY(x,y)) = — (pointwise in x). (25) 

De^oo x 
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Fig. 6. The left panel shows the cross-section of the "UCM extension" profile at 
y = 1.0 for a range of Deborah numbers, De = 0.30, 0.55, 1.00, 1.81, 3.31, 6.03. The 
right panel shows a log-log plot of the width of the peak, both as a FWHM measure 
and for a FENE-P-like definition (see text). 

It is then evident that for given L 2 , the width will approach a constant value 
of 2/L 2 for De -> oo. 



6 Discussion 



One may ask the question to what extent the analysis above is relevant for 
"real" flow profiles. We already argued above that although the details may 
(and will) change, the qualitative behaviour (localization of the extensional 
stress along the separatrix and non-analyticity) will persist. 

The procedure that was performed to obtain the stress profile can be repeated 
for the "exact" velocity field in Refs. [T7j and [T8] . using numerical integra- 
tion where necessary. That would still yield a solution that does not obey the 
momentum balance equation (|2j), but more importantly, it would not quali- 
tatively change the singularity that we find. Mathematically, this happens as 
follows. In the region near the stagnation point x = y = 0, the velocity field 
is expected to have, in lowest order, always terms like in (jHJ). Hence, in this 
regime, the equation always has "zero mode" solutions f(xy 2 ) of the convec- 
tive part of the constitutive equation. Only in special cases will this term not 
be excited. Furthermore, close to the wall, the dominant behaviour of the y- 
independent homogeneous solution must go as exp[— l/(Dey)) and this then 
has to be compensated by a similar essential singularity in y/xy. Thus, close 
to the wall this essentially singular behaviour is robust. Further away from 
the wall, however, the behaviour will change if the velocity field crosses over 
to a different y-dependence. 

Similar considerations apply to other modifications that are used to give a 
more realistic description within the UCM framework. As long as the veloc- 
ity field is reasonably well-behaved, the extension profile on the separatrix is 
determined only by the velocity on the separatrix. 

In fact, we see that the only ingredients that we need in order to obtain our 
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results, are the general form of the velocity field Eq. (J4]), the "advective" part 
of the constitutive equation, as in Eq. ([TBI) , and the absence of a stress diffu- 
sion term in the constitutive equation. This strongly suggests that both the 
localization and the nonanalytic behaviour at x = depend on these ingre- 
dients and not so much on the particular form of the constitutive equations. 
One would then expect that every viscoelastic flow with a stagnation point 
gives rise to singular behaviour downstream of the stagnation point. 

We already mentioned that the stress profile we obtained does not satisfy 
the Navier-Stokes equation (momentum balance). It would be interesting to 
see how we might modify the velocity profile to have the system obey the 
Navier-Stokes equation. It is difficult to do this analytically, but given the 
singular structure of the stress field, it is to be expected that a correction to 
the velocity field will show a similar singularity for intermediate flow rates. As 
we mentioned in the introduction, for sufficiently large flow rates, this coupling 
may cause a breakdown of the steady state flow solution [T2IT3] . 



7 Conclusions 

Not only purely extensional flow but also (reverse) wall stagnation flow shows 
strong localization of extensional stresses and non-analytic behaviour. At mod- 
erate and high Deborah numers, a birefringent strand is formed, with an es- 
sential singularity at its centre. The polymer extension becomes appreciable at 
a finite distance from the wall, and the width of the strand either decreases or 
becomes constant with increasing Deborah number, depending on the defini- 
tion. Asymptotic scaling relations can be found for both the distance from the 
wall and the strand width. Although concrete numerical and analytical results 
have been calculated only for a fixed velocity field, the qualitative behaviour 
that we find should carry over to more realistic flows, at least if the usual 
constitutive equations are to be taken seriously. If an artificial stress diffusion 
term is added to the constitutive equation, as is sometimes done in numerical 
simulations to enhance stability, the mathematical singularity will disappear. 
Nevertheless, the behaviour we have analysed here should still dominate most 
of the profiles if the diffusion constant takes on any value that is physically 
reasonable. 

We suggest that the localization and the singular behaviour are caused by the 
purely advective character (without diffusion) of the constitutive equation, 
combined with the extensional character of the flow. We further suggest that 
the occurrence of this behaviour does not depend on the exact formulation of 
the constitutive equations, as long as these include stress advection and stress 
relaxation. The non-analyticity at the separatrix is a result of the absence 
of diffusion, which essentially decouples the domains on either side of the 
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separatrix. 



As we noted in the introduction, internal stagnation point flows were recently 
shown to exhibit singular behaviour of the type \yf, where y is the upstream 
distance from the singular line. Such singularities may easily cause numerical 
problems since derivatives of sufficiently high order diverge. In the case of wall 
stagnation flow, the situation may at first sight not be as bad for numerical 
approaches, since all derivatives remain finite. Nevertheless, for large Deborah 
numbers there is a rapid crossover region, the width of which decreases as 
1 / De 2 for large De. This region will have to be resolved numerically!!] 
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